function R = compute_residual_image_succ(image, k)
if nargin < 2
    k = 1;
end
R = image;
for i = 1:k
    R = succ_predict_image(R);
end
R = image(1+k:end-k, 1+k:end-k) - R;
end

function p_mean = succ_predict_image(image)
p_rb = predict_image(image, 1, 1);
p_lb = predict_image(image, -1, 1);
p_rt = predict_image(image, 1, -1);
p_lt = predict_image(image, -1, -1);
p_mean = (p_rb + p_lb + p_rt + p_lt)/4;
end

function x = predict_image(image, xdiff, ydiff)
% xdiff, ydiff should be  +1 / -1
[M,N] = size(image);
cA = image(2+ydiff:M-1+ydiff,2:N-1);
cB = image(2:M-1,2+xdiff:N-1+xdiff);
cC = image(2+ydiff:M-1+ydiff,2+xdiff:N-1+xdiff);

maxV = max(cA, cB);
minV = min(cA, cB);

x = cA+cB-cC;
ind = (cC <= minV);
x(ind) = maxV(ind);
ind = (cC >= maxV);
x(ind) = minV(ind);
end